While related to the normal form model, this one fully exapands the higher order terms, thus showing interesting behaviors not present in the Hopf normal form model (see E.W. Large, F.V. Almonte and M.J. Velasco A canonical model for gradient frequency neural networks, Physica D 239 (2010) 905-911):
$$ \dot{z} = z\left(\alpha+j\omega+(\beta_1 +j \delta_1)|z|^2 + \frac{(\beta_2 +j \delta_2)\varepsilon |z|^4}{1-\varepsilon|z|^2}\right) + \frac{x}{1-\sqrt{\varepsilon}x}\cdot\frac{1}{1-\sqrt{\varepsilon}\bar{z}} $$with $x$ being a complex-valued input (could be external or connections to other oscilators in a network of oscilators).
$\mathcal{P}(\varepsilon, x) = \frac{x}{1-\sqrt{\varepsilon}x}$ is a passive non-linearity, and $\mathcal{A(\varepsilon, \bar{z})} = \frac{1}{1-\sqrt{\varepsilon}\bar{z}}$ is an active non-linearity. We can express it as follows:
$$ \dot{z} = z\left(\alpha+j\omega+(\beta_1 +j \delta_1)|z|^2 + \frac{(\beta_2 +j \delta_2)\varepsilon |z|^4}{1-\varepsilon|z|^2}\right) + \mathcal{cP}(\varepsilon,x)\mathcal{A}(\varepsilon,\bar{z}) $$Lets express the external input $x$ and the oscillator state $z$ in polar form:
$$ x=Fe^{j\theta}\\ \bar{z} = re^{-j\phi} $$Then,
$$ \begin{align} \dot{r} &= r\left(\alpha + \beta_1 r^2 + \frac{\beta_2\varepsilon r^4}{1+\varepsilon r^2}\right) + \frac{F(rF\varepsilon+\cos(\theta-\phi)-F\sqrt{\varepsilon}\cos\phi-r\sqrt{\varepsilon}\cos\theta)}{(1-2F\sqrt{\varepsilon}\cos\theta+F^2\varepsilon)(1-2r\sqrt{\varepsilon}\sin\phi+r^2\varepsilon)}\\ \dot{\phi} &= \omega + \delta_1 r^2 + \frac{\delta_2\varepsilon r^4}{1+\varepsilon r^2} + \frac{F(\sin(\theta-\phi)+F\sqrt{\varepsilon}\sin\phi-r\sqrt{\varepsilon}\sin\theta)}{r(1-2F\sqrt{\varepsilon}\cos\theta+F^2\varepsilon)(1-2r\sqrt{\varepsilon}\sin\phi+r^2\varepsilon)}\\ \end{align} $$Is important to note that $x$ and $z$ must satisfy the following conditions:
$$ \begin{align} |x| &< \sqrt{1/\varepsilon} \\ |z| &< \sqrt{1/\varepsilon} \\ \end{align} $$Assuming $\mathcal{P}(\varepsilon, x)\mathcal{A}(\varepsilon, \bar{z} = 0$), we can rewrite the system as:
$$ \begin{align} \dot{r} &= r\left(\alpha + \beta_1 r^2 + \frac{\beta_2\varepsilon r^4}{1+\varepsilon r^2}\right)\\ \dot{\phi} &= \omega + \delta_1 r^2 + \frac{\delta_2\varepsilon r^4}{1+\varepsilon r^2} \\ \end{align} $$
In [17]:
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from pygrfnn import Zparam
from pygrfnn.oscillator import zdot
def rdot(r, z, force):
"""
Compute the derivate of an oscillator's amplitude w.r.t. time, as a function
of the oscillator's amplitude
Args:
r (`array_like`): oscillator's amplitude's to analyze
z (:class:`Zparam`): oscillator parameters
force (`float`): external force
"""
return z.alpha*r + z.beta1*r**3 + z.beta2*z.epsilon*r**5./(1-z.epsilon*r**2) + \
force
z_params = Zparam(0.3, 1.0, -1.0, 0, 0, 1.0) # oscillator params
x = 0 # external force
r = np.arange(0, 1.0, 0.001)
drdt = rdot(r, z_params, x)
plt.figure(figsize=(8,8))
plt.plot(drdt, r)
plt.xlim([-1, 1])
plt.grid(True)
plt.show()
In [ ]: